// atmos_profile.h

#ifndef DODO_ATMOS_PROFILE_H
#define DODO_ATMOS_PROFILE_H

#include <string>
#include <cmath>
#include <stdexcept>
#include <boost/date_time/posix_time/posix_time.hpp>

#include <quigon/basic/constants.h>
#include <quigon/basic/basic_utils.h>
#include <quigon/basic/basic_profile.h>
#include <quigon/basic/basic_serie.h>
#include <quigon/basic/basic_para.h>
#include <quigon/basic/serie_profile_algo.h>

namespace quigon
{


typedef Profile<double, double> Atmos_Profile;
// In a Atmos_Profile, the bin_size & start_value are height(with correction of station height),
// While in a lidar Line, they are distance from lidar.

// 2 versions of load_atmos_profile, for making a profile that matches lidar line.
Atmos_Profile load_atmos_profile(
        const std::string& filename,
        size_t size,
        double bin_size,
        double start_value
);

//requires field in para:
// size_t "size"
// double "bin_size"
// double "start_value"
// double "station_height"
// Link_Stuff<Time_Serie> "elev_angle"
template<class TPara>
Atmos_Profile load_atmos_profile(
        const std::string& filename,
        const TPara& para,
        const boost::posix_time::ptime& t = boost::posix_time::not_a_date_time
)
{
	using namespace std;

	try
	{
		double sin_elev_angle;
		if (t.is_special())
			sin_elev_angle = sin(M_PI / 180.0 * para.elev_angle()()); //wow, elev_angle is a Link_Stuff<Time_Serie>
		else
			sin_elev_angle = sin(M_PI / 180.0 * para.elev_angle()[t]);

		return load_atmos_profile(
		               filename,
		               para.size,
		               para.bin_size * sin_elev_angle,
		               para.start_value * sin_elev_angle + para.station_height
		       );
	}
	catch (...)
	{
		func_err("load_atmos_profile", filename);
		throw;
	}
}

// Notice the Atmos_Profile version of serie_to_profile() is
// result = serie_to_profile<Atmos_Profile> (serie, size, bin_size, start_value);

// cal_beta_m() is calculating backscatter coef of atmos molecules.
// unit is km^-1 / sr
Atmos_Profile cal_beta_m(
        const Atmos_Profile& t_profile,
        const Atmos_Profile& p_profile,
        double lambda = 523.0, // laser wavelength in um
        double rho_n = 0.0284  // need recheck of this value.
) ;

// generate a new profile, that each val is the accumulation of all old_val before.
// e.g. old: 1 2 3 4 5 6,  new: 1, 3, 6, 10, 15, 21
// if with ratio = 2, then new: 2, 6, 12,20, 30, 42
// this is originally designed for accumulating beta_m.
// ratio is a factor that will be multiplied to the result.
// offset indicates that the value is at the center of a bin.

// required Tprof property:
// operator[] for read/write,
// .size(),
// copy-able
template <class Tprof>
Tprof gen_acc_profile(const Tprof& profile, size_t start_index, double ratio = 0.0, bool offset = true)
{
	try
	{
		Tprof result(profile);
		if (offset)
		{
			result[start_index] = profile[start_index] * 0.5;
			for (size_t i = start_index + 1; i < profile.size(); ++i)
			{
				result[i] = result[i - 1] + (profile[i - 1] + profile[i]) * 0.5;
			}
		}
		else
		{
			for (size_t i = start_index + 1; i < profile.size(); ++i)
			{
				result[i] += result[i - 1];
			}
		}
		for (size_t i = start_index; i < profile.size(); ++i)
		{
			result[i] *= ratio;
		}
		return result;
	}
	catch (...)
	{
		func_err("gen_acc_profile");
		throw;
	}
}

template <class Tprof>
Tprof gen_acc_priori_profile(const Tprof& profile, size_t priori_index, size_t start_index, size_t end_index, double ratio = 1.0)
{
	try
	{
		Tprof result(profile);
		if (start_index > priori_index)
			start_index = priori_index;
		if (end_index < priori_index)
			end_index = priori_index;

		result[priori_index] = 0.0;
		for (size_t i = priori_index + 1; i <= end_index; ++i)
		{
			result[i] = result[i - 1] + (profile[i - 1] + profile[i]) * 0.5;
		}
		for (size_t i = priori_index - 1; i >= start_index; --i)
		{
			result[i] = result[i + 1] - (profile[i + 1] + profile[i]) * 0.5;
		}

		for (size_t i = start_index; i <= end_index; ++i)
		{
			result[i] *= ratio;
		}

		return result;
	}
	catch (...)
	{
		func_err("gen_acc_profile");
		throw;
	}
}
}

#endif

